Correlations in optical phonon spectra of complex solids 
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Spectral correlations in the optical phonon spectrum of a solid with a complex unit cell are analysed 
using the Wigner-Dyson statistical approach. Despite the fact that all force constants are real, we 
find that the statistics are predominantly of the GUE type depending on the location within the 
Brillouin zone of a crystal and the unit cell symmetry. Analytic and numerical results for the 
crossover from GOE to GUE statistics are presented. 



Wigner-Dyson statistical analysis pi has become a 
widespread approach to characterise quantum spectra of 
complex dynamical systems, such as nuclei, atoms and 
molecules, disordered quantum electron systems and elec- 
tromagnetic resonances in chaotic cavities [y-gf . Both ex- 
perimental and theoretical studies of such systems point 
to the fact that the quantum energy levels universally 
obey the same mutual correlations and spectral rigidity 
as the eigenvalues of random Gaussian-distributed ma- 
trices. In the present paper, we apply such a statisti- 
cal approach to characterise vibrational spectra of multi- 
component solids with a complex unit cell structure, in 
order to assess the correlation properties of their optical 
phonon spectra taken at various positions in the Brillouin 
zone (BZ). 

In random matrix theory (RMT), different ensembles 
of matrices, which are defined by their fundamental sym- 
metry, result in distinct level correlation properties. The 
two most common in condensed matter physics are the 
Gaussian ensemble of real symmetric matrices (GOE) 
which adequately describes the spectral correlations in 
quantum chaotic electron systems with time reversal 
symmetry, and the Gaussian ensemble of Hermitian ma- 
trices (GUE) which describes chaotic electron billiards 
with time-reversal symmetry broken by a magnetic field. 
One may expect a set of 3A^ coupled classical oscillators 
to be an example of a system with time-reversal sym- 
metry, thus obeying the GOE spectral statistics. This is 
indeed the case with the acoustic spectroscopy of irregu- 
larly shaped solid resonators [g|,0 , or with the spectrum 
of a regularly shaped system consisting of coupled oscil- 
lators whose masses are random [g[. Below we demon- 
strate the counter-intuitive result that, on the whole, the 
spectral correlations of the optical phonon modes associ- 
ated with the same (albeit arbitrary) point of the BZ of 
a complex solid obey the GUE statistics. This complies 
with an earlier observation |Q] on the electronic struc- 
ture of highly excited bands in solids. In particular, we 
report a study of the distribution function P{s) of the 
nearest-level-spacing in vibrational spectra of a complex 
solid based on numerical simulations of crystalline struc- 
tures, both with and without mirror reflection symmetry 
in the unit cell. We also analyse the dependence of sta- 



tistical properties on the phonon wave number Q within 
the BZ, and obtain a detailed description of the crossover 
between the limiting regimes of GOE-type correlations 
specific to the phonon frequencies exactly at the center 
of the BZ where Q = 0, and of GUE-type for sufficiently 
large values of Q. 

To simulate a complex solid, we have adopted the fol- 
lowing model. The unit cell of a crystal was taken in the 
form of a parallelepiped consisting of A^ = 8 x 10 x 12 
atoms with equal pair interactions but randomly cho- 
sen masses arranged on an fee lattice. The unit cell 
size L determines the periodicity of the Bravais lattice 
of the entire solid, L = niLi -t- n2l'2 + n^L^. The spec- 
trum of optical phonons in it can be found by deter- 
mining all 3iV normal modes corresponding to the lin- 
earised set of equations for coupled harmonic oscillators, 
TOjit? = — ^;l i^Tjj (L) uf_,_L, where i, j are atomic posi- 
tions within the unit cell and a, (3 denote cartesian com- 
ponents. Fourier transforming we obtain for each given 
point in the BZ (wave number Q), 



^if(Q) = Efc(^£+ijCij 



'"^^ - 4<5L+ij<5"^), (1) 



where k is the interatomic force constant, /? , , , = 
( j — i — L)"/|j — i — L|, and (j — i — L) belongs to the 
first coordination sphere of an fee lattice. Note that the 
latter implies non-zero matrix elements only for nearest 
neighbors and applies restrictions to the sum over L in 
Eq. (|l|). The sets of phonon frequencies {uJkiQ)} (which 
result in a spaghetti of 3A^ dispersion curves, each corre- 
sponding to a particular phonon branch) can be obtained 
by solving numerically the eigenvalue equation 



det(i:i -uj^ 1)^0 , D{Q) = M-^^^KiQ)M- 



1/2 



(2) 



This equation also defines the dynamical matrix, £'(Q), 
in which the complexity in the composition of a solid 
is introduced by the random diagonal matrix M^J^ — 
niiS^^Sij. For the masses TOj we have used a box distri- 
bution over the interval of masses [(m) — 6m, (m) -|- Sm] 
with 5m/ {m) = 0.3. This corresponds to an r.m.s. of 
mass disorder {Sm.'^y^'^ /{m) ~ 0.17. The mean value (m) 
defines the cut-off frequency Wc = \/8fc/ (to) of the vibra- 



tional spectrum in the disorder-free limit and sets a char- 
acteristic scale to measure the eigenfrequencies, whereas 
lengths are measured in units of the lattice constant. As 
explained below, such a model also allows us to exploit 
an analogy between the numerical results obtained by us 
and the properties of spectra of systems exhibiting quan- 
tum diffusion. To avoid complications brought into the 
problem by the localisation effects related to the vibra- 
tions of light atoms in a heavy matrix, we restrict our 
statistical analysis to the frequency range cj/wc < 0.95. 

At the middle of the BZ (Q = 0) and at other sym- 
metry points, such as Q = (0, 7r/L2,7r/L3), the dynam- 
ical matrix D is real and symmetric. For an arbitrary 
Q 7^ 0, the matrix elements of I?(Q) related to the sites 
on the edges of the unit cell (we use the nearest-neighbor 
interaction) acquire complex phase factors which make 
the whole dynamical matrix complex and Hermitian, 
Z?(Q) = Ds{Cl) + iDA{Ct). where i?s(Q), Da{CI) arc 
real symmetric and antisymmetric matrices, respectively. 
Without imposing any spatial symmetries onto the unit 
cell, these two scenarios span the two limiting cases for 
the D-matrix symmetry; the difference between these two 
regimes manifests itself in the form of the normalized 
distribution function P{s) of the nearest-level-spacing, 
s = {ujk+i — cijfe)/A, where A is the mean level spacing. 




FIG. 1. Nearest-level-spacing distribution of optical 
phonon spectra for a non-symmetric unit cell and various val- 
ues of Q: (a) Q = 0, (b) Q = [0, 7r/4L2, 37r/4L3], and (c) 
Q = [0.05/Li,0,0] (inset : Q = [Q,-n /L2,tt/ L-j.]). Fits to 
histograms are explained in the text. 

We have constructed P{s) employing a two-step aver- 
aging procedure. The first step is to use an ensemble of 
50 random realizations of the distribution of masses in 
the sample. A further averaging of P{s) is applied over 
a broad frequency range, namely, 0.45 < lo/oJc < 0.95, 
for each of the calculated spectra, by using an observa- 
tion that they become statistically indistinguishable af- 



ter having been rescaled by A g. Note that we do not 
analyze the lowest part of the spectrum because of poor 
statistics. The nearest-level-spacing distribution function 
for the optical phonon spectrum at Q = is shown in Fig. 
H(a). It coincides with the random matrix theory predic- 
tion for real symmetric matrices given by Wigner-Dyson 
distribution function for the GOE m. In contrast, opti- 
cal phonons with Q 7^ exhibit a nearest-level-spacing 
statistics which is best fitted by the Wigner-Dyson distri- 
bution function for the GUE. A typical P(s)-histogram is 
plotted in Fig. Jl|(b) , and is compared with the GUE an- 
alytical result m. A similar observation has been made 
earlier by Mucciolo et al |^ in relation to the electronic 
band structures in crystals. For completeness, the inset 
in Fig. ^c) shows the P(s)-histogram for a corner of the 
BZ, which is of GOE-type. Fig. |l|(c) also illustrates the 
form of P{s) in the intermediate regime between two dis- 
tinct statistical classes. It is compared with the result of 
a fit based on the RMT prediction Q for an interpolating 
ensemble between the GOE and GUE which contains a 
single fitting parameter. 

The inter-ensemble crossover takes place at relatively 
small values of Q (the parameter responsible for the 
imaginary part of the dynamical matrix) , as it happens 
with a similar crossover in the energy spectra of chaotic 
electronic systems in the presence of a weak magnetic 
field U^. Note also that for small Q's, the function P{s) 
has the following asymptotic behavior: it resembles the 
GUE distribution function at s ^ 0, whereas it follows 
the GOE analytical result for large s. This is because 
the splitting of levels with s < 1 is more sensitive to a 
small antisymmetric addition to a symmetric dynamical 
matrix than the splitting of rare pairs with s ;3> 1. The 
above fact also implies that the crossover P{s) has a form 
which cannot be simply reduced to a trivial mixing be- 
tween two typical GOE and GUE distributions ||lO| , [ll| . 
The crossover studies were based upon the analysis of 
spectra of 600 random realizations of atomic masses in 
the unit cell for various values of Q along Li. The re- 
sult of numerical analysis of the GOE-GUE crossover in 
different intervals of vibrational spectra is shown at the 
bottom of Fig. g in the form of a gray scale shading of 
the parametric plane of g = Q • Lj^ = QLi and frequency 
uj, where a darker colour stresses higher similarity of P(s) 
to that of the orthogonal symmetry class, and white in- 
dicates the dominance of the GUE spectral statistics. 

Below we present semiclassical arguments which pro- 
vide the crossover value of the rescaled Q, q — QL, as 
function of the frequency of the mode considered. Our 
approach consists of viewing the dynamics as that of 
a wave-packet of lattice vibrations (in our case optical 
modes), spreading over a unit cell treated here as a re- 
gion of disordered medium. This treatment is analogous 
to the semiclassical treatment of parametric spectral cor- 
relations developed in the studies of quantum disordered 
electron systems subjected to a weak magnetic flux [pl|2|. 



In the present analysis, it is assumed that the participa- 
tion of each atom in a given optical phonon mode can be 
described semiclassically for time intervals shorter than 
tn ~ 1/^j after which the discreteness of the spectrum 
for each value of Q starts to dominate. The spread of 
vibrations over the unit cell and the role of the individ- 
ual atoms in this dynamics is considered as a diffusion 
of waves through an fee lattice with mass disorder and 
is determined by the interference pattern of a variety of 
equally probable diffusion paths U^. These paths are 
independent of the exact value of Q, whereas the phases 
of diffusive waves involved in such an analysis are large 
and random, so that one obtains correlated spectra [y,y . 
However, the cyclic (albeit non-periodic) boundary 
conditions impose a torus geometry in the unit cell. If 
W is the total number of windings around the torus then 
the type of correlations depends on the phase factor e^^"^ 
determined by the phase difference between a path with 
a positive W and its time-reversed counterpart with -W. 
Therefore, Sip gives us a measure of time- reversal symme- 
try breaking and is controlled by the 'external' parameter 
q ^ 1. A relevant crossover parameter can then be found 
by estimating the r.m.s. value of the symmetry break- 
ing phase Sep, acquired by waves whose propagation is 
followed along a diffusive path with the maximal length 
allowed by the limit set by the time, tn ~ 1/A [p|,pl|, 
where A r^ 1/L^v and v{u}) is the acoustic phonon den- 
sity of states. 
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FIG. 2. GOE-GUE crossover. The gray-scale plot illus- 
trates the level of deviations of the calculated P(.s) from the 



GOE result, a = 



J[p(.,)-PcOE(s)l^da 



(0 < a < 1). Black 



is used for < cr < 0.25, whereas, the width of the boxes is 
defined by the spectral window used to built P{a). 

The parameter q affects the phase of the partial am- 
plitudes only for paths that cross the unit cell, i.e., paths 
whose 'winding' number is non-zero. Hence, the total 
symmetry breaking phase is S^p ~ qW ^ where W = 'Y^Wi 
is a winding number made of approximately tujtu ran- 
dom contributions Wi = ±1. The time t]j ^ l? jD with 



D the diffusion coefficient, is typical of a diffusive spread 
over the unit cell of a vibration with frequency w. For 
our model, D = ^ (s^) t where (s^(w)) is the angular 
and polarization average of the squared phonon group 
velocity in the clean limit and t~^ is the scattering rate, 
caused by the variation of atomic masses. The latter 
can be estimated |1J] as t^^{lu) ^ uj'^v ( (5'm?'^ / (m) J . 
It follows from our considerations that an estimation of 
the r.m.s. value of 5(p yields {5ip'^)^^^ ^ q^JtiiJtD- The 
crossover can be assigned to {Sip'^y^^ r-^ 1. Therefore, the 
form of the crossover line is determined by 



gc^ VWi^^ VAWF^occ^/y^i^Hy, (3) 

which is shown in the top of Fig. 2 and is in agreement 
with the numerically obtained gray scale plot. 

The correlation properties of spectra of a solid are sen- 
sitive to the geometrical symmetries g of its unit cell 
structure, such as the existence of a mirror plane in it. 
Below, we extend the numerical analysis to a solid with 
n < 3 mirror symmetry planes in the unit cell, each 
characterised by unit vector %. Numerical simulations 
similar to the ones described above (with statistics col- 
lected from 100 random realizations of masses and av- 
eraging extended over the spectral window [0.45,0.95]) 
show that the effect of geometrical symmetries on the 
spectral statistics depends on the point in the BZ. This 
is because the phonon momentum allows for breaking 
both the orthogonal and point-group symmetries. Typi- 
cal forms of the function P{s) for each symmetry case are 
shown in Fig. S, for various wave vectors Q. The solid 
lines in these figures illustrate the RMT result for over- 
lapping sequences of GOE or GUE spectra with equal 
fractional densities 0]. When Q = 0, the spectra split 
into sequences of levels corresponding to different pari- 
ties, so that P{s) coincides with what one would expect 
for 2" overlapping GOE's (e.g., P(s ^ 0) = 1 - 2""). 

For finite Q, spectral statistics is determined by the 
orientation of Q with respect to the mirror symmetry 
planes, falling into one of the GOE or GUE classes. A 
summary of all distinct statistical limits is given in Ta- 
ble 1. Note that for Q with all non-zero components 
perpendicular to a symmetry plane, correlations are of 
GOE-type. This effect is due to invariance of -D(Q) un- 
der the combination of reflection and complex conjuga- 
tion operations, which results in a real representation of 
the dynamical matrix B. Thus, for a unit cell with 3 
mirror symmetry planes only the orthogonal ensemble 
statistics is realised for an arbitrary Q. The crossover 
parameter qc which determines which class from Table 1 
should be expected, can be estimated using Eq. (||). 
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FIG. 3. Nearest-level-spacing distribution of optical 
phonon spectra for various unit cell symmetries and values 
of Q: Q = (B, E, L), Q = [7r/4Li,0,0] (A, G, M), 
Q = [0,0,7r/4L3] (D, H), Q = [7r/4Li, 37r/4L2,0] (C, F, P), 
Q = [0,37r/4L2,7r/4i3l (K),andQ = [37r/4Li,7r/3L2,7r/4i3] 
(J, N). Histograms A, F, N(C, J) are compared to the 
GOE(GUE) analytic result. Other fits are the P(s) for a 
number of overlapping RMT spectra with equal fractional 
densities: two GOE(GUE) in B, G, P(D, K), four GOE(GUE) 
in E, M(H), and eight GOE(the Poisson distribution is also 
plotted) in L. 
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TABLE I. Summary of all distinct statistical limits of op- 
tical phonon spectra for n < 3 mirror symmetry planes in the 
unit cell, each characterised by unit vector fji. RMT ensem- 
bles preceded by a number denote a sequence of corresponding 
overlapping spectra with equal fractional densities. Each box 
contains the condition for these statistics to be realised at 
some point Q in the BZ. 
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